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Abstract — We investigate the density-gradient (DG) transport 
model for efficient multi-dimensional simulation of quantum con- 
finement effects in advanced MOS devices. The formulation of the 
DG model is described as a quantum correction to the classical 
drift-diffusion model. Quantum confinement effects are shown to 
be significant in sub-lOOnm MOSFETs. In thin-oxide MOS capac- 
itors, quantum effects may reduce gate capacitance by 25% or 
more. As a result, the inclusion of quantum effects in simulations 
dramatically improves the match between C-V simulations and 
measurements for oxide thickness down to 2 nm. Significant 
quantum corrections also occur in the I-V characteristics of short- 
channel (30 to 100 nm) n-MOSFETs, with current drive reduced 
by up to 70%. This effect is shown to result from reduced inver- 
sion charge due to quantum confinement of electrons in the chan- 
nel. Also, subthreshold slope is degraded by 15 to 20 mV/decade 
with the inclusion of quantum effects via the density-gradient 
model, and short channel effects (in particular, drain-induced 
barrier lowering) are noticeably increased. 

I. Introduction 

I N THE fast-moving electronics industry, the classical drift- 
diffusion (DD) model of electron transport has steadfastly 
remained the dominant model for industrial numerical simula- 
tion of electronic devices [1], 50 years after the model’s first 
description [2], and 35 years after Gummel [3] described a 
robust numerical solution method. This uncommon longevity 
stems partly from momentum (due to familiarity and experi- 
ence), but it derives mainly from the fact that, with tuning for a 
given technology, the DD model continues to provide a combi- 
nation of computational speed, robustness, and acceptable 
accuracy which has been difficult to match with other models. 
However, there is a growing realization that technologists (not 
just researchers) can not ignore quantum effects much longer. 
Any refinement or replacement of the DD model should ideally 
maintain all of its benefits while including quantum effects. We 
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describe such a model and initial promising simulation results 
in this paper. 

Clearly, the semiconductor industry is in a period of feverish 
advancement, with new generations of electronics technology 
being developed every 2 years, rather than every 3 years as pre- 
dicted by recent history [4]. For the industry-dominant MOS- 
FET, gate lengths and oxide thicknesses of production devices 
will shrink towards 50 nm and 1 nm respectively over the next 
decade [5]. This and other work shows that quantum confine- 
ment effects will significantly affect the operation of such 
ultra-small devices. It is unclear how well additional fitting 
parameters in the classical DD model can account for these 
quantum effects. However, formulating alternative, computa- 
tionally efficient, accurate, and robust multi-dimensional elec- 
tronic device models including quantum effects has been very 
challenging. Two approaches have been taken In the attempt to 
meet this challenge: (i) employing full quantum models for 
conventional device simulation, and (ii) adding quantum cor- 
rections to classical models such as DD. 

The first approach includes the use of quantum models such 
as non-equilibrium Green’s function [6] and Wigner function 
[7]. The Green’s function model contains a high level of quan- 
tum mechanical and scattering detail, but in multi-dimensions 
it is many orders of magnitude more costly than the DD model. 
The Wigner function model trades quantum mechanical and 
scattering detail for somewhat better computational efficiency 
in comparison to the Green’s function model, but it suffers 
from unsolved numerical robustness problems [8], and is still 
much too costly for production-level computations, even in 2- 
D. We note that 1-D simulations can provide some information 
about electronic device operation, but 2-D (or 3-D) simulations 
are essential for sufficient detail and accuracy for industrial 
application. 

The second approach for including quantum effects in device 
simulations is to add quantum corrections to classical models. 
For the DD model in a MOSFET, approximate quantum cor- 
rections have been applied to the inversion charge profile [9], 
the silicon band gap near the oxide [10], and both the band gap 
and the gate oxide thickness [11]. Alternatively, a 1-D 
Schrodinger computation has been added to the DD model to 
account for quantum confinement effects in the inversion layer 
[12], [13]. These approaches may be computationally efficient, 
but they assume a particular device structure (MOSFET with 
planar gate) and operating region (inversion). In contrast, quan- 
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turn corrections to the hydrodynamic model [14], [15] are gen- 
eral in terms of allowed device structures and operating 
regions, but they suffer from significant disadvantages in com- 
putational robustness and cost versus the DD model [16]. 

In this work, we describe a device simulation model which 
uses the second approach to meet the multi-dimensional quan- 
tum challenge without the disadvantages of the approaches 
described above. This model is the density-gradient (DG) 
quantum correction to the DD model [17]. The DG model adds 
quantum confinement and (optionally) tunneling to the DD 
model in a general, compact, and computationally efficient 
manner. This yields a model which meets all of the require- 
ments of a replacement for the DD model. In this work, we 
focus on quantum confinement effects in MOS devices, using 
an implementation of the DG model which does not include 
quantum tunneling. After describing the model, the remainder 
of this paper compares classical and DG model simulation 
results. First, we compare capacitance- voltage (C-V) curves for 
thin-oxide MOS capacitors, showing a significant reduction in 
capacitance due to quantum confinement of electrons and 
holes. We then show the large reduction in drive current of an 
ultra-small MOSFETs when quantum confinement is included. 
We also show the degradation of small MOSFET subthreshold 
slope and the increase in short-channel effects (DIBL). We 
conclude with a summary of these results. 


n. Density-Gradient Model 
The drift-diffusion and density-gradient models of carrier 
transport in an electronic device can be written identically, 
including Poisson’s equation and the electron and hole trans- 
port equations: 


V-(eV\|/) = -p = -q(p - n + C); (la) 

§7 = ^ = V + DVn) ; (lb) 

= V(pn p Vv p + D F Vp)-, (lc) 


where \j / is electrostatic potential, n and p are electron and 
hole densities, C is fixed charge density, £ is permittivity, p is 
total charge density, q is electron charge, J is current density, 
and [i and D are mobility and diffusivity of the respective car- 
riers. 

In the classical DD model, the electron and hole “drift poten- 
tials” are just the electrostatic potential: =\|/ . In the DG 

model, y n and \\f p have quantum corrections: 
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The expressions for the quantum potentials \\f qn and \\f qp are 
derived from the Schrodinger equation, based on the finite cur- 
vature (energy) and strict continuity of wavefunctions [17], 
[18]. Since the quantum potentials represent a net effect for all 
wavefunctions, they do not incorporate quantum mechanics 
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Figure 1: Basic n-MOSFET structure and biasing. Bias V G between gate 

and source contacts controls current / D flowing between source and drain. 
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Figure 2: Schematic of classical (dashed lines) and quantum (solid lines) 

carrier density profiles of a vertical cut through an n-MOSFET in inversion. 
Classical densities are sharply peaked and discontinuous at the silicon/oxide 
interface, and can not penetrate into the oxide. Quantum profiles are smoothly 
peaked below the interface, continuous at all interfaces, and can penetrate (tun- 
nel) into the oxide. 

exactly into the DG model, so r n and r p may be used as fitting 
parameters. In this work, we take r n =r^=3, which is the high 
temperature limit [19], [20]. The quantum potentials act to 
smooth the carrier density profiles by reducing their second 
derivatives (curvature). In fact, the DG model forbids disconti- 
nuities in the carrier density profiles. 

We focus on quantum confinement effects in the industry- 
dominant device, the silicon MOSFET (Figure 1). The effects 
of quantum smoothing on carrier profiles in a MOSFET operat- 
ing in inversion are depicted in Figure 2. Classical carrier den- 
sities (dashed curves) change abruptly at the oxide interfaces 
from some large external value to zero in the oxide. Quantum 
mechanical carrier densities (solid curves) can not change 
abruptly. The densities must be continuous across the oxide 
interfaces, resulting in significant differences in classical and 
quantum carrier profiles near the interfaces, and penetration of 
quantum carrier densities (quantum tunneling) into the oxide. 

Returning to the DG model in (1) and (2), the five solution 
variables are \|/, n, p, , each of which involve second- 

order partial differential equations (PDEs). The quantum 
potentials \|f and \| / accomplish the expected carrier pro- 

file smoothing, and so will be largest near the oxide interfaces 
where the classical density discontinuities will be smoothed 
out. In this work, we ignore carrier densities (and thus tunnel- 
ing) in the oxide, and we set the carrier densities n and p to 0 
in this region. Thus, only the source-less Poisson equation (la) 
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is solved in the oxide, while the full DG model is solved in 
both the gate and substrate. In order for the carrier densities to 
be continuous across the oxide interfaces, they must approach 
zero just outside the oxide. With this and the usual boundary 
conditions (BCs) for the DD model, all of the BCs of the DG 
model are well defined except for those for the quantum poten- 
tials at the oxide interfaces. 

The two quantum potentials have relatively large, unknown 
values at the oxide interfaces, since this is where they act to 
force electron and hole concentrations smoothly to zero. [Actu- 
ally, the quantum potentials would have to be infinite at the 
oxide interfaces to force the carrier densities to exactly zero 
there, so we instead set the densities at the interfaces to a small 
but non- zero value; I0’ 8 /cm 3 was used in this work.] Thus, a 
suitable constraint on the values of the quantum potentials at 
the oxide interfaces (a Dirichlet BC) is not available. By the 
same reasoning, enforcing a Neumann BC on the (unknown) 
gradients of the quantum potentials is also not valid. One solu- 
tion to this lack of quantum potential BCs at the oxide inter- 
faces is to solve the entire five-PDE model in the oxide as well 
as in the adjoining silicon and poly gate. This also implements 
tunneling in the DG model. In this case, the discontinuity in the 
quantum potentials is determined by the silicon-oxide band 
offsets, and the gradients would be continuous across the inter- 
face. A model implementing this approach is being developed 
[ 21 ]. 

Another solution to the boundary condition challenge is to 
use the quasi-Fermi (QF) model [22] of carrier transport, 
which can be described as a change of variables from the DD 
model. In the QF model, the continuity equations are: 

fj = (3a) 

j, - ■ V <W, T V < 3b) 

where, including DG quantum corrections, the QF energies 
are: 

0„ = \\f-(kT/q)\n(n/n i )^y\f gn ; (4a) 

0 p = V+ikT/q^nip/nJ + y^, (4b) 

and n i is the intrinsic carrier concentration of the semiconduc- 
tor. In (4a) and (4b), we have assumed a Maxwell-Boltzmann 
energy distribution of the carriers. Analysis using Fermi-Dirac 
statistics, which is more exact but more expensive and difficult 
to implement, will be presented in the future. 

Note from (3a) and (3b) that at the interface between a semi- 
conductor and an insulator, the electron and hole QF energies 

6 and 0 in the semiconductor have zero gradient normal to 

~n p v 

the interface, since current flow J n , J into the insulator is 
zero. [Recall that tunneling current is zero in this work.] The 
DG version of the QF model can use these constraints on 0^ 
and 0 as the additional BCs needed to solve the quantum 
potential PDEs. Before writing the final DG model used for 
this work, we note that only steady-state simulations were 
used, so the time derivatives were eliminated from the continu- 
ity equations. Finally, the model solved in this work is: 


(eV\| i) + q(p - n + C) = 0, 

(5a) 

V • (nn„V0 n ) = 0, 

(5b) 

-V • (pn p V0 p ) = o. 

(5c) 

J>'V qn - 2b n W2 J' 1 = °’ 

(5d) 

Jpv q p + 2b p vl 'Tp = 0- 

(5e) 


The five solution variables for these PDEs, in order, are 
\jj, 0 p , n , p . Expressions for the quantum constants b n and 

b p are given in (2a) and (2b). From (4a) and (4b): 

Vgn = 0 n + (^7’ / fl)ln(«/n J .)-V; (6a) 

V qp = $ p -( kT/ <l) ]n (p/«,-)-V- (6b) 

Note that the DG model in (5a) - (5e) is generic, in the sense 
that it can be applied to any electronic device structure. By 
contrast, other quantum-DD models often incorporate struc- 
ture-specific and localized quantum corrections [9]-[ 13], and 
may not allow tunneling to be included. 

For this work, we used fixed mobilities of = 1500 cm W s 
and }X p = 500 cm 2 /Vs in (5b) and (5c), which are roughly 
equal to the intrinsic values at room temperature for silicon 
[23]. [Exact current predictions are not being pursued here; 
only quantum effects on device behavior. More accurate mobil- 
ity models are under development.] For the quantum constants 
in (2a) and (2b), we used m n = QA9m 0 (light electron mass in 
silicon) and m p = 0A9m o (heavy hole mass in silicon). These 
values result in a good match between DG simulations and 
experiment over a wide range of MOS device structures (see 
Section III. A). All simulations were for devices at room tem- 
perature (300K). 

We implemented the above density-gradient model in a par- 
tial differential equation (PDE) solver called PROPHET [24]. 
This simulation tool provided several advantages over the tra- 
ditional approach of coding a numerical simulator specifically 
for the DG model. It allowed us to specify the transport model 
in a script file at a high level to a general-purpose (but highly 
efficient) PDE solver. Also, PROPHET has the necessary facil- 
ities to solve the model in 1-D, 2-D, or 3-D for any specified 
device and test regime. Thus, the PDE- solver approach allows 
for the rapid investigation of a wide range of device structures, 
transport models and physical effects. Indeed, three signifi- 
cantly different versions of the DG model were investigated for 
this work. Using the conventional, model-specific approach, 
programming all of the major and minor model variations 
investigated would have taken many times as long. The classi- 
cal QF model was also implemented in PROPHET, and was 
used for all classical model simulations. 

III. Results and Discussion 

A. Thin Oxide MOS Capacitors 

The switching efficiency of a MOSFET is largely determined 
by its gate capacitance, which measures the ability of gate 
biases to control the carrier density and thus current flow below 
the gate oxide (see Figure 1). Thus, it is critical for simulations 
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Position (nm) 

Figure 3: Classical and DG quantum-corrected carrier density profiles in a 

MOS capacitor operating in inversion. The inset details the electron profiles on 
a linear scale, showing the dramatic quantum correction at the oxide interface. 



Figure 4: MOS Capacitor low frequency C-V curve comparison for (100 

pm)^ area, 2. 1 nm thick oxide. The DG model reproduces measurements (dots; 
courtesy of H-P Labs) much more closely than the classical model. 

to accurately predict gate capacitance. The quantum repulsion 
of carriers from both gate and substrate oxide interfaces, as 
depicted in Figure 2, makes the oxide appear to be typically 1- 
2 nm thicker than it is. This effect is already quite noticeable in 
state-of-the-art commercial products, which have gate oxide 
thicknesses as low as 2 nm. The effect will quickly increase 
over the next decade, with gate oxide thicknesses predicted to 
shrink to 1 nm or less by 2012 [5]. 

To test this prediction, 1-D MOS capacitors with oxide 
thicknesses from 2 to 8 nm were simulated, and C-V curves 
were compared to those from classical simulations and from 
experimental measurements of the same structure [25]. As an 
example, the electron and hole densities for a 2.1 nm oxide 
device biased in inversion are shown in Figure 3. The resulting 
carrier density profiles were as predicted in Figure 2: near zero 
at the oxide interfaces, with the inversion or accumulation 
charge peak 0.5 to 1.5 nm beneath the Si-oxide interface, rather 
than exactly at the interface as in the classical model. Figure 4 
compares C-V curves (capacitance versus gate bias V G ) for 
the 2.1 nm oxide MOS capacitor. As expected, the DG model 



Oxide Thickness (nm) 


Figure 5: Percent error in simulated capacitance versus gate oxide thickness 

at V G = -2V (accumulation). The DG model maintains accuracy at least down 
to 2. 1 nm. The accuracy of the classical model deteriorates rapidly for oxide 
thicknesses below 4 nm. 

reproduces measured data much more accurately than the clas- 
sical model for this very thin oxide. 

To summarize the results of the C-V simulations over the 
range of oxide thickness considered, Figure 5 shows the frac- 
tional error in simulated capacitance (compared to measured 
data) versus oxide thickness for the classical and DG models. 
To simplify the plot, a single gate bias of V G = -2V (accumu- 
lation) was chosen, since this condition is most critically 
affected by quantum effects and is least affected by other 
unknown parameters such as the poly doping level [26]. Here 
we see that the DG model maintains accuracy at least down to 
2 nm, while the accuracy of the classical model deteriorates 
rapidly for oxide thicknesses below 4 nm. 

B . Short Channel MOSFET 

As discussed in Section I, 1-D simulations can provide only 
limited knowledge of device operation. 2-D and 3-D simula- 
tions are often required for an accurate analysis of the opera- 
tion of state-of-the-art (highly non-planar) devices. For 
example, the MOSFET operating current I D , which flows hor- 
izontally from source to drain, is largely controlled by the ver- 
tical electric field from the gate. Thus, modeling operating 
current accurately requires at least a 2-D analysis. Our previous 
work [26] was the first to show that 2-D simulations are quite 
feasible using the DG model. The fact that the DG model is 
also general (not structure-specific) allows it to work without 
modification or tuning for complex (e.g., non-planar) struc- 
tures. 

Figure 6 shows the simulated drain characteristic (drain cur- 
rent versus drain bias at a series of gate biases) for a very 
aggressively scaled MOSFET with a 30 nm gate length and 2 
nm gate oxide thickness. This device approximates the state of 
the art (or slightly beyond) in research labs [27]. Computed 
drain characteristics for both the quantum-corrected DG model 
and the classical quasi-Fermi model are shown. At each gate 
bias, the DG current is 20% to 70% below that predicted by the 
classical model. This represents a serious decrease in the cur- 
rent drive capability of the device due to quantum effects. We 
want to point out, however, that even this aggressively scaled 
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Figure 6: Simulated drain characteristics for 30 nm gate length, 20A gate 

oxide MOSFET. Results for classical and DG simulations are shown for 0.4V 
< Vq < IV. The classical current is 25-60% larger. 



Figure 8: Subthreshold characteristic of 30-nm MOSFET. Quantum effects 

degrade the OFF state at low bias, the ON state at operating biases, and the sat- 
uration effect at high bias. The net effect is like a parallel conductance. 



Figure 7: Relative drain current and channel charge (at mid channel with 

V 0 s = IV) for the DG model (DG result divided by classical result). The result 
indicates that, in the normal operating range, the DG current reduction is dom- 
inated by reduced channel charge, with only minor quantum transport effects. 



Figure 9: DIBL (drain-induced barrier lowering, or threshold voltage reduc- 

tion) versus channel length. The channel of the 30-nm MOSFET was stretched 
to create longer channel MOSFETs. DIBL is somewhat more severe (stronger 
threshold variation with channel length) when quantum effects are included. 


30 nm MOSFET still performs like a switch. Figure 6 shows a 
negligible simulated drain current at logic-zero gate bias (0 V), 
and significant drain current at logic-one bias (1 V). Admit- 
tedly, the device is in need of extensive engineering to mini- 
mize short-channel effects. 

One question which needs to be answered is whether the 
reduced DG current in Figure 6 is due to the reduction in chan- 
nel charge, quantum transport effects along the channel, or 
both. Figure 7 compares the channel charge and drain current 
from the DG simulation (relative to the classical values) at full 
drain bias (IV) for the 30 nm MOSFET simulated in Figure 6. 
The close match between these curves over the full range of 
gate biases indicates that the DG current reduction is domi- 
nated by reduced channel charge, with only minor quantum 
transport effects. It makes sense that quantum effects are minor 
in the transport direction in the DG model, since this model 
only significantly affects the potential and carrier profiles near 
abrupt heterojunctions and insulating interfaces. The relatively 
smooth potential in the transport direction results in small 
quantum potentials in this direction, and correspondingly small 


quantum effects on current. We note that other quantum mod- 
els which include the effect of discrete quantum energy levels 
in the channel may predict more significant quantum effects 
due to carrier transport along the channel. 

As final confirmation that the DG model and this tiny 30 nm 
MOSFET work correctly, Figures 8 and 9 compare classical 
and DG simulations of the subthreshold characteristic and 
DIBL (drain-induced barrier lowering versus channel length) 
for this device. In Figure 8, the subthreshold slope is 90 to 92 
mV/decade of current for the classical model, and 105 to 1 10 
mV/decade for the DG model. Thus, quantum confinement 
effects significantly degrade the subthreshold slope. 

In Figure 9, DIBL was simulated by determining the thresh- 
old voltage versus channel length at maximum drain bias ( V D 
= IV). To create longer channel MOSFETs, the 30 nm device 
was stretched at the center of the channel. The threshold volt- 
ages was computed as the gate bias at which the potential bar- 
rier to current flow between source and drain was the negative 
of the built-in potential in the substrate (0.494V in this case). In 
the classical model, this barrier is the classical potential at the 
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gate oxide interface. In the DG model, the classical and quan- 
tum potentials combine to produce the barrier between source 
and drain, and the minimum barrier to carrier flow is slightly 
beneath the oxide surface. Figure 9 shows that the DG model 
shifts the threshold voltage by about 150 to 200 mV. Further, 
DIBL is somewhat worse (threshold voltage varies more rap- 
idly) with the DG model. In both cases, variation with channel 
length is very severe near 30 nm for this device. Thus, although 
the MOSFET works, inevitable process variations would likely 
make this device structure unsuitable for ULSI. 

IV. Computational Efficiency 

The forgoing simulation results show that the DG model 
makes it feasible to include quantum effects accurately and 
generally in multi-dimensional electronic device simulation. In 
this section, we go further to show that the DG model is in fact 
quite efficient in accomplishing this. In particular, we compare 
the computational cost of the classical and DG models for the 
simulations in Section III. 

For the MOS capacitor simulations of Section III. A, we used 
a non-uniform 1-D grid with 320 to 450 points (depending on 
oxide thickness) for both classical and DG simulations. On an 
SGI 02 workstation, a typical C-V curve trace (81 bias points) 
took 40 seconds for the classical model, and 74 seconds for the 
DG model. For this case, computation time including quantum 
effects with the DG model is less than a factor of 2 larger than 
that for purely classical simulations. 

For simulating MOSFET I-V curves in Section III.B, we 
used an identical 2-D grid with about 1750 points for both the 
classical and DG simulations. On the same workstation, a typi- 
cal I-V curve trace (51 bias points) took 452 seconds with the 
classical model and 2383 seconds for the DG model. We note 
that the DG model in Section II proved to be unstable at times, 
in which case we used a slightly modified model with Jn and 
Jp as solution variables rather than n and p . An unresolved 
error in the Jacobian for the modified DG model results in lin- 
ear convergence (rather than quadratic), and a correspondingly 
longer simulation time. Thus, computing an I-V curve with the 
modified DG model requires typically 5400 seconds. 

These DG model computation times should be compared to 
the orders of magnitude increase in computation time for more 
rigorous quantum models, such as those based on the 
Schrodinger equation [28] or Green’s functions [29]. Since the 
DG model is only moderately more computationally demand- 
ing than the associated classical models, it can even be feasibly 
solved in 3-D [30]. More importantly, the DG model leverages 
all of the tuning and optimization of the industry standard, 
classical drift-diffusion model. Thus, the DG model provides 
practical insight into quantum effects in ultra-small electronic 
devices without the uncertain accuracy or meticulous tuning 
effort that face more rigorous quantum models. 

v. Conclusions 

In summary, we presented the density-gradient as a compu- 
tationally efficient means of including quantum effects in 
multi-dimensional electronic device simulation suitable for 
future MOSFET technology. Computation time is typically 2 


to 10 times that for the purely classical model, making it possi- 
ble to run 2-D electronic device simulations routinely on a 
workstation. With the simulation of numerous 1-D and 2-D 
MOS devices, we also demonstrated the robustness of the DG 
model. In MOS capacitor simulations, classical model predic- 
tions rapidly diverge from measured results for oxide thick- 
nesses below 4 nm, while the DG model maintains error below 
a few percent down to 2 nm. 

In simulations of a 30 nm gate length MOSFET, quantum 
effects are predicted to reduce current drive by up to 70%. 
According to the DG model, this current reduction is almost 
entirely due to the reduced inversion charge resulting from ver- 
tical quantum confinement, while horizontal quantum transport 
effects along the channel are minimal. Finally, we showed that 
quantum effects degrade the subthreshold slope of this small 
MOSFET by 15-20 mV/decade, and increase DIBL. We con- 
clude that the inclusion of quantum effects is essential for the 
accurate simulation of 30 nm scale MOSFETs, and that the 
density gradient model is an efficient way to accomplish that. 
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